The field theory of symmetrical layered electrolytic 
systems and the thermal Casimir effect 

DS Deanf and RR Horgan$ 

f Laboratoire de Physique Theorique, CNRS UMR5152, IRSAMC, Universite Paul 
Sabatier, 118 route de Narbonne, 31062 Toulouse Cedex 04, France 
$ DAMTP, CMS, University of Cambridge, Cambridge, CB3 OWA, UK 

Abstract. We present a general extension of a field-theoretic approach developed in 
earlier papers to the calculation of the free energy of symmetrically layered electrolytic 
systems which is based on the Sine-Gordon field theory for the Coulomb gas. The 
method is to construct the partition function in terms of the Feynman evolution kernel 
in the Euclidean time variable associated with the coordinate normal to the surfaces 
defining the layered structure. The theory is applicable to cylindrical systems and its 
development is motivated by the possibility that a static van der Waals or thermal 
Casimir force could provide an attractive force stabilising a dielectric tube formed 
from a lipid bilayer, an example of which are t-tubules occurring in certain muscle 
cells. In this context, we apply the theory to the calculation of the thermal Casimir 
effect for a dielectric tube of radius R and thickness S formed from such a membrane 
in water. In a grand canonical approach we find that the leading contribution to the 
Casimir energy behaves like —ksTLKc/ R which gives rise to an attractive force which 
tends to contract the tube radius. We find that kc ~ 0.3 for the case of typical lipid 
membrane t-tubules. We conclude that except in the case of a very soft membrane 
this force is insufficient to stabilise such tubes against the bending stress which tend 
to increase the radius. We briefly discuss the role of lipid membrane reservoir implicit 
in the approach and whether its nature in biological systems may possibly lead to a 
stabilising mechanism for such lipid tubes. 



PACS numbers: 87.16.Dg, 05.20.-y 
1. Introduction 

In an recent short communication we reported on a calculation which investigated the 
the possibility that a static van der Waals or thermal Casimir force could provide an 
attractive force across a tube formed from a lipid bilayer, so leading to its stabilisation. 
In this paper we give the details of the general theory of symmetrically layered 
electrolytic systems which underlies that calculation, and explain the details of the 
calculation applying the theory to cylindrical geometry and to a model for the lipid 
bilayer tube. Whilst the motivation for developing the theory presented below is the 
analysis of the the Casimir force in the context of a dielectric tube immersed in water, the 
theory is applicable to any sufficiently symmetrical system consisting of layer containing 
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electrolyte. The Coulomb properties of such systems are described by a Sine-Gordon 
field theory and a full analysis in the case of flat layers has been done with the approach 
which is generalised in this paper. In particular, it allows for the perturbation series 
for the thermal Casimir force to be developed in terms of the dimensionless coupling 
g = Ib/Id where and Id are the Bjerrum and Debye lengths, respectively. 

The behaviour of systems composed of layer of varying dielectric constants was 
first studied by Lifshitz and coworkers [1] and has been subsequently revisited by a 
number of authors [2, 3, 4]. The formalism developed is an elegant way of taking 
into account van der Waals forces in a continuum theory. Two types of van der 
Waals forces are accounted for in theses theories, firstly zero frequency van der Waals 
forces whose nature is purely classical and secondly the frequency dependent ones 
due to temporal dipole fluctuations. In terms of thermal field theory the former 
correspond to the zero frequency Matsubara frequency and the latter to the non zero 
frequencies. In order to calculate these latter terms we require information about the 
frequency dependence of the dielectric constants, where as the former only requires 
the static dielectric permittivity. The quantum Casimir effect corresponding to the 
modification of the ground state energy of the electromagnetic field has been intensively 
studied in the case of idealised boundary conditions in a variety of geometries including 
spheres and cylinders [5]. The thermal Casimir effect investigated here has a similar 
mathematical structure though the corresponding effective spatial dimension is one less 
in the calculations. The temperature dependence of the full Casimir in a simplified model 
of a solid dielectric cylinder (and sphere) has been recently examined using a heat kernel 
coefficient expansion [6] . In our analysis of the diffuse limits we make use of summation 
theorems for Bessel functions which were introduced for the study of the Casimir energy 
for cylinders with light- velocity conserving boundary conditions [7]. In this paper we 
will calculate only the zero frequency contribution, also known as the thermal Casimir 
effect. The thermal Casimir effect may also be calculated in the presence of electrolyte 
and the technique we develop here for electrolytic systems within the Debye Hiickel 
approximation is valid in the domain of weak electrolyte concentrations. There is an 
extensive literature on the thermal Casimir effect for systems of layered geometries, 
both without added electrolyte and within the Debye Hiickel approximation [2, 3, 4]. 
Recently the calculation of the thermal Casimir force for layered films at the first order 
of perturbation theory about the Debye Hiickel theory was carried out [8], suggesting 
the possibility of strong non perturbative effects. 

In section 2 we discuss the model for the lipid bilayer tube and review the outcome 
of the calculation applied to this model; in section 3 we present the general theory for 
calculating the free energy of a general symmetrically layered electrolytic system; in 
section 5 we apply the general theory to dielectric layers with cylindrical geometry; in 
section 6 we present the calculation of the Casimir force for the particular case of a 
dielectric tube of thickness 5 and radius R immersed in water; in section 7 we evaluate 
the Casimir force for physically reasonable values of 5, R and in section 8 we present 
some conclusions. 
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2. The lipid tubule 

The behaviour of lipid bilayers is of crucial importance in biophysics. Lipid bilayers 
in water exhibit a huge variety of geometries and structures and in the context 
of cell biology even more varied structures are exhibited. In order to understand 
where biological mechanisms such as molecular motors and cytoskeletal structures are 
determinant in the stability of biological structures, one must first understand the role 
of the basic physical interactions in systems that contain only lipid bilayers, i.e. model 
membrane systems. There has been much study of lipid bilayer shape and elasticity 
using standard continuum mechanics [9, 10]. This basic approach is also complemented 
by more microscopic studies based on lipid structure and lipid-lipid interaction models, 
this approach is of course ultimately necessary to fully understand the physics of bilayers. 
The bilayer is composed of two layers of lipid each layer having the hydrophilic lipid 
head at the surface where it is in contact with water, the interior is composed of the 
lipid's hydrocarbon tails. This layer geometry is stable due to the hydrophobic nature 
of the hydrocarbon tails. Given this non-homogeneous structure one can immediately 
see that a simple continuum elastic sheet type model may have difficulty in predicting 
the mechanical properties of bilayers. 

In certain muscle cells, structures known as t-tubules are found. These are basically 
cylindrical tubes whose surface is composed of a lipid bilayer. Similar structures may 
also be mechanically drawn off bilayer vesicles. The stability of these tubular structures 
requires an explanation. The basic continuum theory[9, 10] predicts that the free energy 
of a tube of length L and radius R is 

F B (L,R) = k -S^, (1) 

where the above expression is strictly speaking the excess free energy with respect to a 
flat membrane of the same area A = 2-kRL and the subscript B refers to mechanical 
bending. Various experimental and theoretical estimates for Kb can be found in the 
literature [10, 11] and they lie between 3 and 30. Note that our definition of Kb differs 
from that used traditionally in the literature k' b by a factor of ir, Kb = ttk,' b . The 
values of k b depend of course on the composition of the bilayer and on the experimental 
protocol used to measure it. One crucial element in both theoretical and experimental 
determinations of Kb is whether the tube is attached to a reservoir of lipid or not, i.e. 
whether the statistical ensemble is grand canonical or grand canonical. Clearly if there 
is no reservoir then any increase in the surface area of the tube will lead to a less dense 
lipid surface concentration, in this case water may be able to become in contact with 
the internal layer composed of the hydrophobic heads and thus a significant increase 
in free energy. If upon, changing the area of the tube, lipids can flow into the tube 
to maintain the local optimal packing then the free energy cost will be substantially 
different. This bending free energy is positive and hence the preferred thermodynamic 
state is the flat one. Mechanical models for membranes vary in their predictions for the 
dependence of Kb on the membrane thickness 5. The models most compatible with the 
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available experimental data predict that 

k b = K A 5 2 /a (2) 

where a depends on the precise model and is generically 0(10) and Ka is the area 
compression modulus [10, 12]. Experimental fits of k b with respect to the membrane 
thickness are compatible with 

K B = K A (5-5 ) 2 /a (3) 

where Sq ~ lnm is an offset necessary to fit the data. We note that when lipid tubes 
are drawn from a vesicle the mechanically applied tension can of course overcome this 
free energy barrier. A natural question motivated by the fact we see these structures 
in cells is whether there are any other mechanisms that could lead to their formation 
and explain their stability. A possible explanation is that electrostatic effects involving 
surface charges and ions (salt) in the surrounding medium could play a role[13, 14, 15] 
. Certain experiments [14] however revealed a relative insensitivity of some system to 
the concentration of salt. There are however other systems where the salt concentration 
does appear important in determining the stability of the tubules [16], in these systems 
the lipid head groups are highly charged. Another explanation has been put forward in 
terms of the geometry of the lipid, notably the tail having a structure such that there 
is a preferred orientation of the tails next to each other, giving rise to a chirality which 
allows the stabilisation of the tubes [17, 18, 19, 20, 21].. This explanation would however 
depend on a more or less mono-disperse lipid bilayer in order to permit this liquid crystal 
like phase. Cell membranes are composed of a wide variety of lipid types and addition 
have proteins present and so it is possible that another mechanism is responsible for the 
stability of these structures 

We adopt a continuum model where the lipid bilayer is modelled as a layer of 
thickness 5 ~ 5 — 10 (nm) and of dielectric constant e' ~ 2eo- The surrounding water 
is also treated as a dielectric continuum of dielectric constant e 80e . We shall also 
adopt a model where the lipid tube is fixed at each end to a flat lipid reservoir and thus 
work in the grand canonical ensemble; this is shown schematically in Fig. 1. 

In this paper we find that the thermal Casimir effect gives a contribution to the 
excess free energy above the flat plane of 

F (L,fl) = -*5^=e (4) 



with 



2 



A 

KC = 64 



3 log ( — J +61og2 + 3 7£ -4 



+ A 4 5(A) , (5) 



where a is a microscopic cut-off corresponding to the molecular/ lipid size below which 
the continuum picture of the dielectric medium breaks down. We note that the sign of 
F c has exactly the same functional form as the bending free energy F B but is of opposite 
sign, meaning that this force tends to collapse the tube and thus helps to stabilise the 
system against the bending energy. We shall show later that with reasonable physical 
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parameters kq ~ 0.5 — 1.0. Thus the Casimir attraction is not able to overcome the 
repulsion due bending as it is predicted by current theories and data. However this 
result is important for several reasons: 

• We show that the Casimir attraction tends to stabilise the tube structure. 

• The presence of the microscopic cut off in k c shows that the physics is ultimately 
dominated by the short scale or ultra-violet physics. This means that weak 
electrolyte concentrations will have little effect on the system as seen in experiments 
given that there are no strong surface charges. 

• We see that Fq and Fb have the same functional form at large R and that the 
behaviour of kc is regulated by the microscopic physics. This means that our 
calculation can be interpreted as a renormalisation of Kb due to the thermal Casimir 
effect. 

• Further attractive,or tube stabilising, interactions may be generated by the high 
frequency Matsubara modes. 

3. The Schrodinger kernel for separable systems 

The mathematical tool that allows us to derive the free energy for electrolytic systems 
with symmetrical layered films is the functional Schrodinger kernel which evolves the 
Sine-Gordon scalar field from some in initial surface to a final surface. Our method 
is applicable when the Laplacian is separable in the natural coordinates describing the 
surfaces bounding the layers of the system. In the surface between the bounding layers 
the electrostatic and chemical properties of the system are uniform i.e. the dielectric 
constants and electrolyte concentrations are constant. It is in this sense that we describe 
such a system as symmetrical. In this case, in D-dimensions, the coordinates can be 
denoted by (x, a) where x is a list of (D — l) coordinates for surfaces a = constant. The 
i-th surface of an N- layer system is described by a = <Ji, where the Oi are constants with 
<7j + i > <7j, < % < N with (Tq and being respectively the minimum and maximum 
values in the range of a. The local electrochemical properties of the system thus depend 
solely on the coordinate a. Our example in this paper will be that of coaxial cylinders 
in D = 3, where x = (9, 0) and a = r, the radius. However, the theory is more general 
than for cylindrical or spherical coordinates, and so we lay the theory out below in 
a general notation but refer to the cylindrical case for clarity where appropriate. The 
dynamics of the field 0(x, a) are defined by its evolution in the Euclidean time coordinate 
t, — oo < t < oo which is given in terms of a. The volume measure is dv = J(a)d<jdx 
and the Euclidean time t(a) is defined by 



For example, in the cylindrical geometry a — r, t — log a and in the planar case 
t — a — z. 




(6) 
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To derive the general form for the kernel it is convenient to express the contribution 
from one layer of the system to the total free energy in dimensionless variables. In a 
previous paper [8] we derived an expression for the grand partition function of a layered 
system in a dimensionless form, and in the present context the effective action is the 
Sine-Gordon field theory which defines the kernel and is written as 

1 / * + M 

'Via Aug Jv 12 

where the region defining the layer is bounded by two neighbouring surfaces Si and S2 
defined, respectively, by a — a±, a = a 2 , (J\ < cr 2 , and has volume Vi2- All lengths are 



S = ~f dv (V0) 2 + ^ I dv cos( v ^0), (7) 
8n JV12 4:7iq Jv 12 



measured in terms of the Debye length, Id = 1/m, where m = ^2pe 2 (3/e is the Debye 
mass, p is the ion density of the bulk reservoir to which the electrolyte solution within the 
layer is connected, and e is the dielectric constant in the layer. The other fundamental 
length in the theory is the Bjerrum length, lg = e 2 (3/4ire, and the dimensionless coupling 
constant is given by g = Ib/Id- The dimensionful field is given in terms of by the 
rescaling 

— ^ ■ (8) 

The renormalisation constant Z(g) is associated with the ion chemical potential p 
conjugate to p, and removes the divergences due to the unphysical charge self- 
interactions. In Eq. (7), p has been substituted by the reservoir density p using the 
relation 

" = Z[a)p ■ Z(g) = ( ' (9) 

where the above subscript B indicates that expectation value is for an infinite bulk 
system. 

The total partition function is constructed as a convolution of the kernels of the 
layers in sequence, and to carry this out the dimensionful description must be restored. 
For multiple layers the action is a sum of similar terms each associated with a layer of 
the system bounded by an inner and an outer surface. In particular, the innermost and 
outermost surfaces are at er and a N corresponding to t — —00 and t — 00, respectively. 
It was shown in [4, 22] that for planar interfaces the Schrodinger kernels which are 
bounded by one or the other of these surfaces are given in terms of the ground-state 
wave-function of the appropriate free Hamiltonian, and that this is sufficient to ensure 
that the overall charge neutrality constraint is respected. In the more general case, 
where the interfaces are non-planar (cylindrical, for example), the Hamiltonian depends 
explicitly on the Euclidean time t and so there is no interpretation in terms of stationary 
eigenstates. However, in the limit t — > ±00 the relevant kernels are separable in the 
boundary fields, and this leads to the same result. 

The action S in Eq. (7) can be decomposed as 

S = + S (0) + AS , (10) 

Aug 
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where V is the volume of the layer and the first term is the ideal contribution. The 
term is the action for a free or Gaussian field theory and is given by 

SV> = ~ I dv [(V0) 2 + 2 1 . (11) 
The interacting part of S is expressed as a perturbation 



AS = — ( dv 

Aug Jv 



Z(g) ( cos(^0) - 1) + 9 4- 



(12) 



and the action Sb for the equivalent bulk system is given by 

S B = --L / dv (V0) 2 + ^ / dv cos( v ^0) , (13) 
8n jb Ang Jb 

which may be decomposed in the same manner as for S. 
The Schrodinger kernel for the layer is defined by 

£(<fe(x),<7 2 ;<Mx),<7i) = t 'VtPe 3 ^, (14) 

where <fo(x) = 0(x, af), % — 1,2, are the boundary values of the field 0(x, cr) on the 
bounding surfaces 5j, respectively. 

In this section we concentrate on the calculation of JO o )(0 2 (x), a 2 ; 0i(x), af) defined 

by 

^(^(x),^;^),^) = f^V<Pe s(0) ^. (15) 

J<t>\ 

The explicit evaluation of for the specified geometry gives the Casimir-effect 

contribution from the layer to the free energy Q = — A^TlogH, where S is the grand 
partition function for the system, and forms the basis for a perturbative expansion of 
Q in terms of the interaction coupling strength g. For an iV-layer system the grand 
partition function for the free theory, E,(°\ is given by the convolution over layers as 

r N 

i=0 

where i(o"o) = —oo, t(a^) = oo, and where the are re-expressed in terms of 

the original, dimensionful, boundary fields so that their values match correctly on the 
common interface separating successive layers. The Casimir free energy is then given 
by 

F c = fi {0) - . (17) 

Here is the equivalent bulk contribution of an independent set of pure bulk systems 
having the same volume and properties as the layers composing the system. In this 
way the generalised force corresponding to the position of any interface is a disjoining 
pressure. 

We shall now show how to explicitly compute i^°)(02, &2\ <f>i, erf) in its dimensionless 
form. The volume measure in Eq. (15) is dv = J(a)dadx. where J (a) is the Jacobian of 
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the measure. Since the functional integral defining is Gaussian in form we explicitly 
find the classical field C which minimises the action by solving the linear field equation 

-(V • J(a)V)fa + J(a)fa = , (18) 

with boundary constraints 

c (x,(7i) = 0i (x), c (x,o- 2 ) = 02 (x) . (19) 

We assume that the operator V • J(<r)V is separable which allows us to write this field 
equation as 



~J(^<Pc - J(a)(V^+l)0 c = 0, 



(20) 



where is self-adjoint and may depend on a but not on derivatives with respect to a. 
The orthonormal eigenfunctions of — V£ are denoted X(s, x) with eigenvalue are A(s, a): 

-V*X(s,x) = A(s,(7)X(s,x) , (21) 

where s is a set of D — 1 quantum numbers. The classical field fa(x, a) is expanded on 
the complete set of functions {X} as 

c (x,a) = ]T T(s,a)X(s,x) , (22) 



where T(s, a) satisfies the ordinary differential equation 

T(s,a) = 0. 



(23) 



We denote two solutions of this equation by i*i(s, a) and F 2 (s, a), where Fi(s, a) is finite 
as t(a) — > — oo and F 2 (s,a) is finite as t(a) — > oo. In addition, these functions with 
different quantum numbers s are orthogonal with respect to the appropriate measure. 
The Wronskian is given by the identity 

J(a)[F 1 (s,a)Fi(s,a)-F{( S ,a)F 2 (s,a)} = 1. 

Then we can write 



T(s,a) = ai(s)Fi(s,a) + a 2 (s)F 2 (s, a) . 
The boundary fields fa on the surfaces Si of the system can be expanded as 
&(x) = ]T Q(s)X(s,x) , <2 < N . 



(24) 
(25) 



(26) 



For the generic layer under discussion we consider the bounding surfaces to be Si and S 2 . 
Comparing with Eqs. (22) and (25), we find the relation between c(s) = (ci(s), c 2 (s)) 
and a(s) = (ai(s) , a 2 (s)) to be 



c = a ■ F(s,a 2 ,a 1 ) , F(s,<r 2 ,<7i) 



/ Fi(s,<ti) 
V ^2(s,ai) 



^1(8,(72) \ 

^2(s,a 2 ) / 



(27) 
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Now using the classical field in Eq. (22) and the definition of from Eq. (11) we find 
that the free classical action is (<f> c ) is given by the boundary term 



S (o) (0c 



8tt J 



J(a)0 c (x,o-)- 



(28) 



where we have used integration by parts. 

We use the expansion of c (x, cr) in Eq. (22) in terms of the coefficients a(s) and 
the expansion of 0i(x), i = 1,2, in Eq. (26) in terms of the coefficients c(s), and also 
use the fact that the functions of the basis set {X(s, x)} are orthonormal. We can then 
eliminate a(s) in favour of c(s), and find from Eq. (28) that 



with 

D = F~ X G , 



G(s,(7 2 ,(7i) 



/ -7(^)^(8,(7!) 

v -J(ai)^(s,ai) 

Then we have 

£(°)(02(x),<72;&(x),<7i) = 1] ^ ( ° ) (S,C 2 (S), ( T 2 ;C 1 (S), ( T 1 ) 



7(^)^(8,(72) \ 

J((7 2 )F 2 '(s,a 2 ) ) 



(29) 



(30) 



K (0) (s,c 2 (s),o- 2 ;ci(s),o-i) = A(s,<7 2 ,<7i) exp (~^ c ( s ) • D(s,a 2 ,cr 1 ) ■ c(s] 



(31) 



where the normalisation factors A(s, cr 2 , o^) arise from the Gaussian integration over the 
fluctuations £(x, t) of the field 0(x, t) about the classical solution. We have that 



Dt, exp I /'"'' do .1(a) 



r\ 2 



(32) 



A(s,(7 2 ,<7i) = 

where 

0(x, o") = c (x, cr) + £(x, a) , (33) 
0(s,cr) = c (s,<7) +|(s,a), 

and where generically we have defined the transform, /(s, cr), of a function /(x, cr) by 
/(s,cr) = |rfx/(x,a)X(s,x). (34) 

The boundary conditions are 

^x,^) =f(x,cr 2 ) = =S> 1(8,(7!) =|(s,cr 2 ) = 0, Vs. (35) 

Then we have 

A(s, a 2 ,a ± ) oc (det[£ CT (s)])- 1/2 , £ a (s) = - ± J(a)-jL + J(<r)(A(s) + 1). (36) 

The determinant can be calculated by diagonalising C t (s) on a basis of orthonormal 
eigenfunctions which satisfy the boundary conditions on £(cr) given in Eq. (35). Whilst 
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yielding the correct result this is not the quickest way to compute A(s, a 2 , erf). The 
Pauli-van-Vleck formula tells us that 



2vr 



det 



d 2 SW(<j) c 



901(902 

Using the expression for S^°\(f) c ) in Eq. (29), we find 



1/2 



A(s,a 2 ,cr 1 ) 



\D 12 (s,a 2 ,(T 1 ) 



2tt 



(37) 



(38) 



(40) 



The Pauli-Van Vleck formula can be derived by analytically continuing the 
Euclidean time variable t to Minkowski time r, by performing the Wick rotation 
t — > t — —it. For this purpose we consider the kernels K and as functions of 
U = t(<7i) rather than <7j, and then K((p 2 , ir 2 ; <f>i, irf) as defined by Eq. (14) is a unitary 
operator which means that 

K^(s,c 2 ,it;c u it)) = Jdc' (k^(s, c',ir'; c 2 , it))* K^(s, c',ir'; c u it) = <5(c 2 - Cl ),(39) 
for any r'. From Eq. (31) 

K (0) (s,c 2 ,iT 2 ;c 1 ,iT 1 ) = A(s, ir 2 , in) exp (-^ c ( s ) • ^ J (s,r 2 ,ri) • c(s)^ , 
where D^s, r 2 , n) = —iD(s,iT 2 ,iTi) is a real symmetric matrix. Then 
\A(s,iT 2 ,iTi)\ 2 J dc' exp Q 2£)( 2 (s, ir 2 , in)(c 2 - ci)c' + £> 22 (s, ?r 2 , in)(c 2 - c?) 

= *(C2-Ci). (41) 

This equation determines A(s, it 2) irf) and, on analytic continuation back to Euclidean 
time t and re-expressing as a function of a, we find that A(s,a 2 ,a±) is given by Eq. 
(38). 

The kernel (<j) 2 (x), iT 2 ; 0i(x), irf) analytically continued to Minkowski time r 
describes the time evolution of the wave-function in the associated quantum mechanics 
problem. The Hamiltonian associated with K^(s, c 2 , t 2 ; c±, ti) is 



W(s,c,f) = -i^ + J 2 (a)(A(s )( x) + l)c 2 , 



(42) 



where a = a(t) is defined by inverting Eq. (6). This Hamiltonian contains an explicit 
time dependence and so the usual quantum mechanical analysis becomes more general. 
The Euclidean version of the Schrodinger equation for wave-function -0(s, c, t) is 

-J^(s,c,t) = H(s,c,t)4>(s,c,t) . (43) 

This equation is also satisfied by K^(s, c, t; d ', t') regarded as a function of (c, i) for fixed 
(c',t'). As remarked earlier, because the Hamiltonian is explicitly t-dependent there are 
no stationary states associated with 7i(s,c, t). However, in the limit that either t — > oo 
or t — > — oo the kernel X^(s, c, i; c', £') is a separable function of c and d except in 
one particular case. These properties will be elucidated in the context of cylindrical 
interfaces discussed in the next section. The connection between the grand partition 
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function in statistical mechanics and the related quantum mechanical formulation is the 
one outlined above between the imaginary and real time formalism ([23]). 

The final outcome for the contribution from quantum numbers s to the kernel for 
the free field theory in the layer, up to an irrelevant factor, is 

fr (0) (s,C2,CT2,ci,o-i) = ?i = 1 == exp ■ D(s, a 2 , af) ■ cj , (44) 



\H(s,a 2 ,cr 1> 
where, using Eqs. (27) and (30), we find 

1 

£>(s,<7 2 ,<7i) = 



/ W(s,<7 2 ,<7i) 



H(s,a 2 ,af 



\ 



\ 1 W(S,<J 1 ,<7 2 ) J 

We have used the identity for the Wronskian in Eq. (24), and have defined 

W(s, (Tj, of) = J(a i )[F 1 (s, <Tj)F 2 (s, a*) - F[(s, ai)F 2 (s, af)\ , 
H(S, <Tj, <7i) = Fi(s, <Ti)F 2 (s, of) - F 2 (s, (Ti)Fi(s, of) . 



(45) 



(46) 



4. Concentric cylinders 

We now apply the formalism of the previous section to the case of two concentric 
cylinders of length L in the z-direction and radii r\ and r 2 , respectively , with r\ < r 2 . 
The separable coordinates are x = (9, z), a = r and the Euclidean time coordinate 
is t — log(r) (note, all coordinates are considered dimensionless at this stage). Hence, 
r — > 0(oo) =>• t — > — oo(oo). In what follows, we work with r rather than t for simplicity. 
The volume measure is 



dv = rdrdOdz 
and Eq. (21) becomes 



J(r) 



1 &■ 



r 2 d6 



rj2 

2 + — )X(s,^,^) = A(s,r)X(s,M, 



(47) 



(48) 



with solution 



y ' ' ; 2tt 



s = (n, p) , n G Z, — oo < p < oo , 
X(s,r) = (n 2 /r 2 +p 2 ) . 

Eq. (23) is then 

A A ^2 



(49) 



d d fit ■ i 
-— r— + — + (j/ + l)r 



T(s,r) 







,2 

7 

This is Bessel's modified equation [24], and the two required solutions are 
Fi(s,r) = /„(Pr), F 2 (s,r) = X n (Pr) . 



(50) 



(51) 
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where P 2 = p 2 + 1. We note that the Wronskian condition of Eq. (24) is satisfied by 
these solutions since 

Pr[I n (Pr)K' n (Pr) - I n {Pr)I' n {Pr)\ = 1 . (52) 

Then, using Eqs. (45) and (46), we find 

/ W n (Pr 2 ,P ri ) 1 \ 



£>(s,r 2 ,ri) 



H n (Pr 2 ,P ri ) 



\ 1 W n (P ri ,Pr 2 ) ) 



(53) 



where 



W^Pr^Pn) = Pr % [I n (Pr ] )K' n {Pr. l ) - I' n (P ri )K n (P rj )] , 

H n (Pr v P n ) = I n (Pn)K n (P rj ) - K n {Pn)I n {Prj) . (54) 

The contribution from quantum numbers s = (n,p) to the kernel for the free field theory 
in the layer, up to an irrelevant constant factor, is then 

K (0 \s,c 2 ,r 2 ,c 1 ,r 1 ) = 1 _ exp (~\c- D(n,p,r 2 ,n) ■ cj . (55) 

''\H n (Pr 2 ,Pr 1 )\ v 1 J 



4-1. Asymptotic behaviour 

We use the definition of K^(s, c', r', c, r) in Eq. (55) and the asymptotic behaviour for 
the Bessel functions given in Eq. (A. 3) to derive the behaviour of as r — > and 
r — > oo. Because we need to consider the case when the Debye mass m is zero we carry 
out the analysis using dimensionful coordinates. This follows easily if we interpret p 
and r as carrying dimension, with P 2 = p 2 + m 2 , and rescale D — > (5eD. 

4-1.1. r' — > 00 In the limit r' — > 00 (t' — > 00) the natural boundary condition for the 
scalar field is = 0, which corresponds to d — 0, and we impose this condition from 
now on. For the various cases we find 

p — rn = : 

/ (3e c 2 \ 



K^\n = 0,0, c' = 0,r';c,r) ~ exp 



^log(r'/r) V 2 log(r'/r) 



1 / /?e c 2 \ 

:ex P — IV • ( 56 ) 



y/F—t™*\ 2(t'~t) 

This is the free particle kernel for Euclidean time ([23]). It is the one case where 
the kernel is not separable in the c, c' variables. The important feature of this result 
for the charge neutrality condition is that 

^log(^ (0) ) = -fcjtbfi - as t' ^ 00. (57) 

Also 

K i0) (n > 0,0, c' = 0,r';c,r) ~ V^i exp (-^ log r'/r^j exp [^^ nc ^j ■ ( 58 ) 
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This is the harmonic oscillator ground-state in c with associated energy E = n/2. 
The prefactor contains the correct (Euclidean) time- dependent factor exp(— n/2(t' — 

*))■ 

P = y/p 2 + m 2 > : 

We define the function V n (z) by 

^> = " ' < 59 > 

Then 

K (0) {n,p,c' = 0,r'; c,r) 

~ (27rPr') 1/4 exp f -ip/) -— L= exp f-^K(Pr)c 2 ^) . (60) 

The related Schrodinger equation, which has a time-dependent Hamiltonian, satisfied 
by K^(s, d, r'; c, r) considered as a function of c and t = log(r), is 

-^(s,c,t) = (-^ + ^(^ ! + Ii V))# ) c li ). (61) 

It can be verified that the different forms listed above in the limit t' — > oo do, indeed, 
satisfy this equation. 

4-1.2. r — > In the limit r — > (i — > — oo), for the various cases, we find 
p = m = : 

1 / /3e(c-c') 2 \ 

This is the free particle kernel for Euclidean time ([23]). As before 

^-log(K^) -> as i -> -oo. (63) 
which we shall show ensures charge neutrality. Also 



K {0) (n > 0,0, c',r';c,r) ~ v 7 ^ exp log r'/r^j exp |^-yn(c 2 + c' 2 ) 



(64) 



This is the harmonic oscillator ground state in both c and d with associated energy 
Eo — n/2. The prefactor contains the correct (Euclidean) time- dependent factor 
exp(-n/2(i' - *)). 

P = VP 2 + rn 2 > : 

We define the function U n (z) by 

U n (z) = ^ . (65) 
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Then 



1 { f3e. 



K<M(n = 0,p,c',r';c,r) ~ exp -^[/ (P r ')c' 2 , (66) 

X (0) (n > 0,p,c',r';c,r) ~ 

exp Qnlog(Pr)^) exp ^-y nc ^ ^ = = ex P (^^U n (Pr')c 2 ^j . (67) 

As before these asymptotic forms satisfy the related Schrodinger equation. 
4-2. The cylindrical membrane 

The i-th layer of a system of N concentric cylindrical layers has volume Vi and is bounded 
by the cylindrical surfaces Si and S^+i, where So is the innermost surface and Sn is the 
outermost. The Debye mass and dielectric constant associated with the bulk reservoir 
connected to the i-th layer are denoted by 6j, respectively. The tube is of length L in 
the z-direction which is parallel to the symmetry axis of the cylinders. The contribution 
to the grand partition function from this system is the convolution 

- N-l 

E M = / Vc J[ J[ ir i (s,Q+i(s),r m ;Q(s),r i ) , (68) 

J s i=o 

where s = (n, p) as before and the field measure is 

N 

Vc = nn dc *( s ) > ( 69 ) 

s i=0 

with boundary condition cjv(s) = V s. 

The grand partition function for the whole system including the bulk reservoirs to 
which the different layers connect is 

^ N-l 

3 = ^ = II Sfl,(Vi) , (70) 

i=0 

where S#. (V^) is the bulk grand partition function for % — th layer of volume V$. This 
can be calculated using the bulk action defined in Eq. (13) with chemical potential \ii 
and dielectric constant for a torus of volume V*. 

The free energy is then F = — /c B TlogS and the forces acting on the interfaces 
and the stability of the system can be deduced from F. From Eq. (11) and the earlier 
discussion the perturbation theory for F can be developed as a loop expansion with 
expansion parameters gi = m^B where rrii is the Debye mass of the i-th. layer and Ib 
is the Bjerrum length. The expansion is a cumulant expansion about the quadratic, 
or free, field theory which is described by the quadratic approximation to the grand 
partition function, S( ) of the system 

-.(0) _ 

^m= /©cnn 1 ^ (0) (s,Q + i(s),r m ;c l (s),r l ), (71) 
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WATER 



Figure 1. Horizontal cross section through idealised tubule configuration shown are 
dimensions and dielectric permittivities 

and with denned in terms of the free- field action as described just above. Each term 
in the product over s on the RHS of Eq. (71) has an exponent which is a quadratic form 
in the interface field variables c(s), where the notation c(s) has been used to signify the 
vector of all the cfs) < i < N associated with the given set of quantum numbers 
s. The integral over the boundary fields with respect to the measure T>c is therefore 
Gaussian, and can be done exactly. The free energy = — ksT log S 1 - ^ contains 
the ideal gas contribution and this one-loop term. The one-loop term consists of a 
contribution from the normalisation factors of the Kf \ and from the determinant of 
the matrix defining the quadratic form in the exponent which arises from the Gaussian 
integration over the boundary field values. The Casimir forces acting on the system 
are determined by the one-loop contribution. In [4] the attractive Casimir forces acting 
between the faces of a planar soap film were discussed and derived in this manner, 
and the contribution to two-loop order in the cumulant expansion of the interaction 
AS defined in Eq. (12) for the planar film was presented in [8]. In general, the loop 
perturbation theory can be carried out in the same way for any symmetrical layered 
electrolytic system such as that constructed from concentric cylinders or spheres. This 
perturbation theory will be pursued in a future publication. In the next section we 
analyse the case of a thin cylindrical membrane for which N — 3. 

5. The Casimir force for a dielectric tube 

In Figure 1 the cross-section of a tube of inner radius R formed from a membrane of 
thickness 5 is shown, with radii for the boundary surfaces defined to be 

5 5 

r = 0, n = R--, r 2 = R+-, r 3 = oo , (72) 

with 5 < R. 

In this paper we concentrate on the Casimir force acting on the tube described 
above and shown in Figure 2 in which the electrolyte densities, and hence the Debye 
masses, are zero in all three layers, the membrane is of fixed thickness S with dielectric 
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E 



- FLAT MEMBRANE 
RESERVOIR 



Figure 2. Vertical cross section through idealised tubule configuration showing tubule 
of length L and radius R bridging two flat bilayers 



constant e 2 = = 2e vac , and the inner and outer layers are filled with water so that 
e o — e 3 — e w — 80e mc . The Casimir force is thus due purely to the discontinuity in the 
dielectric constants at the membrane surfaces and is a function of the radius R of the 
inner cylindrical layer. We shall show that the Casimir force in this case is attractive, 
tending to collapse the tube. The tube is of length L which we assume is large on any 
relevant scale, and that there is a reservoir of membrane so that the tube radius R 
can change without the membrane needing to stretch. For example, the system can be 
thought of as made from a flat sheet of membrane on to which the tube connects and 
which acts as a reservoir of membrane as the tube expands or contracts as shown in 
Figure (2). Alternatively, the membrane can be folded at one end of the tube and act as 
a reservoir. Both scenarios are possible in biological systems where the membrane is a 
lipid bilayer, although for the latter it is difficult to calculate the free energy of a given 
volume of lipid in the reservoir so presents a problem of normalisation. However, it 
should be emphasised that this picture may nevertheless be an important feature of the 
stability of lipid tubules and needs further analysis. In either case, because we assume 
that the membrane is not stretched as the area of the tube increases, there is no elastic 
energy stored in the tube except that due to the curvature and the surface area of the 
system is constant. In what follows we shall assume a reservoir of flat membrane as 
shown in Figure (2). 

As an intermediate step we define the grand partition function, E$ W (R, 5), for the 
membrane normalised by that of an equivalent water filled region, and its associated 
free energy F$ W (R, 5) by 

S£V0M) = ^ D iR,6) - , F$UR,S) = -k B T\og^Sw{R,5), (73) 

where is given by Eq. (71), and where is the grand partition function of a 
system filled with water only: = for 5 = 0. For the grand ensemble with a 
reservoir consisting of a flat membrane of the same thickness, shown in Figure 1, we 
must subtract the free energy of a flat membrane of equivalent area to the tube. We 
then find that the free energy appropriate to calculate the Casimir force due to the tube 
geometry is 

F C (R,5) = F$l(R,6) - 2ttRLF oc (5), (74) 
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where 



lim 



(75) 



2ttRL ' 

where F^S) is the free energy per unit area of flat membrane of thickness 5. Using 
the expression for E^(R,8) in Eq. (71) and the asymptotic expressions for as 
r — > in Eqs. (64) and (65) it can be seen that the dependence on co(s) , which are 
the coefficients determining the boundary field value 0o( x ) at r = 0, cancels out Vs 
between numerator and denominator in Eq. (73). This is independent of whether the 
integrals over the Co(s) are done or not; it is a consequence only of separability in the 
limit r — > or, in the case n — 0, that becomes independent of the c (s). Thus we 
may set co(s) = Vs in what follows and omit the integrals over the co(s). 

The difference between the numerator and denominator in E^ W (R, 5) in Eq. (73 is 
then due to the different contributions of the membrane layer between radii r\ = R — S/2 
and r 2 = R + 5/2 which has dielectric constant €m in the numerator but ew in the 
denominator. We then have 

Z { $ W (R,5) = / Y[d Cl (s)dc2(s) g (0) (s,c 2 (s),r 2 ; ,Ci(s),ri)K 1 (0) (s,c 2 (s),r 2 ; ,ci(s),ri),(76) 

J s 

where 



Q (0) (s,c 2 ,r 2 ; ,ci,ri) 



lim 

r ^0 



if (s, 0, r 3 ; c 2 , r 2 )ir 2 (s, ci, r x ; 0, r ) 



(77) 



^(s,0,r 3 ;0,r ) 

where in the dielectric constant in the denominator kernel is ew From the 
asymptotic expressions in Eqs. (56-67) we find that has a simple form for 

(n,p)^(0,0) (c.f. s = (n,p)) 

Q {0 \s,c 2 ,r 2 ;,c 1 ,r 1 ) 



n l/2 



2ir\K n (pr 2 )I n (pri 



cxp 



-/3e 



K(pr 2 )c 2 exp 



/3e 



w 



U n {wi)c\ 



(78) 



where U n ,V n are defined in Eqs. (66) and (59), respectively. Note, that all Debye 
masses are zero here. For large argument U n , V n — > 1, and so both Gaussian forms 
are convergent and integrable. takes the form of a normalised product of a 

generalisation of harmonic oscillator ground state wave-functions in C\ and c 2 , and 
satisfies the appropriate Schrodinger equations in these variables. For large R and 
n > these functions become the usual ground-state oscillator wave-functions with 
moo = f3e w n, which agrees with the analysis of the planar film of [4]. 
For (n,p) = (0,0) we find 



Q(°)(s = 0,c 2 ,r 2 ; ,ci,ri) = 
where t j = log r\ and iy = ^ 



lim 

?'3 — *oo 
r ^0 



(3e 



w 



t 



30 



1/2 



cxp 



3 ' 



2vr t 32 tio 

Also, we have that 




cxp 



/3e 



w c i 



2 t 



(o). 



0,c 2 (s),r 2 ;,ci(s),ri) = 



/3e 



27T t 21 



1/2 



exp 



(3t w (ci - c 2 ) 5 



<-21 



,(79) 



(80) 
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and so we find the contribution from s = mode to the partition function is 

ll/2 



0,R,5) 



lim 

r_3 — *oo 
r ^0 



GW^Mt30 



1 . (81) 



fwi^W — ^u)t 2 \ + ^W^mH{}_ 

Thus, the zero mode does not contribute to the free energy of the membrane. However, 
it is the relevant mode to show that the charge neutrality condition holds. The total 
charge operator £ is given by 

£ = ew j dxE(x, r 2 + r]) — ew J dxEfari — rj) , (82) 

where -E(x, r) is the radial component of the electric field and rj is a small positive 
length; thus, we measure the field in the water just outside the membrane surfaces. We 
have that 



d 



1 <JUcm)> = 



-<7r(x,t)) ,(83) 



J(r(t)f3e w 

now considering E and as functions of t — log(r). Here 7r(x, t) is the momentum 
operator conjugate to 0(x, t) and is given by 7r(x, t) = /3e0(x, t) using standard theory. 
The Schrodinger representation of 7r(x, t) [22] then gives 

1 , . 5 



£(x,t) = 



J(t)/?e 



(Si) 



j(t)/r^' 1 50(o) 



50(x) 
1 o 

10, t) , 



|0, t). 



The contribution to (E) from the integral over the surface at r = r\ is then 

1 ,~ , 5 



where 0(0) is the zero-mode field. Using Eqs. (57) and (79) we then find 

e w (0,ti| 0(0) |0, if) 



(Si 



lim 

to->- oo J{ti) 







(84) 



(85) 



(86) 



t\ — 1 

A similar result holds for (£ 2 ), the contribution to (E) from the surface integral at 
r = r 2 ; thus we find that (E) = 0. This analysis can be repeated to show that all 
moments of E vanish: (E n ) = 0, Vn > 0; it is this condition that ensures charge 
neutrality of the system. 

We now calculate the free energy by summing over all s = (n,p) mode 
contributions. From the previous section we find 

^i (0) (s,c 2 (s),r 2 ; ,ci(s),ri) 



(3e 



M 



_2-K\H n (pr 2 ,pr 1 )\ 
where as before c = (ci, c 2 ) and 



D(s,r 2 ,n) 



1/2 



cxp 



c(s) • D(s,r 2 ,n) • c(s) 



(87) 



H n (pr 2: pri) 



' W n (pr 2 ,pr 1 ) 1 N 

y 1 W n (pr u pr 2 ) J 



(88) 
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and W n and H n are defined in Eq. (54). Using Eqs. (76), (78), and (87) we find 



/il^ d Cl (s)dc 2 (s) 





1 

2 



where 



and 





1/2 




1/2 


\K n (pr 2 )I„(pr 1 )\_ 




_\H n (pr 2 ,pr 1 )\_ 




z(s) ■ E M (s,r 2 ,rf 


•c( 


a)) , 


(89) 


,r 2 ,n) + e w X(s 


r2,n) , 


(90) 



X(s,r 2 ,ri) = diag[C/ n (pri), K»(pr 2 )] , 
is a 2 x 2 diagonal matrix. A more useful alternative expression is 



(91) 



1/2 / rfc ( s ) ex P c ( s ) ' E m(s, r 2 , n) ■ c(s)^J 
s J dc(s) exp ^-^ c ( s ) • E w (s, r 2 , n) ■ c(s)^J 



,(92) 



where 



E w (s,r 2 ,rf) = e w (D(s,r 2 ,r 1 ) + X(s,r 2 ,r 1 )) . 



(93) 



The denominator is the contribution from a pure water-filled system. The two 
expressions for E^ W (R,5) are the same since using the Wronskian identity, Eq. (24), 
we find 



det[D(s, r 2 , rf) + X(s, r 2 , n)] 
Then we have 



[K n (pr 2 )I n (pr 1 )H n (pr 2 ,pr 1 )] 1 . 



(94) 



B(s,r 2 ,ri) 



1 fe M V' 2 



2 VevK/ 



ndet[B(s,r 2 ,n)] , 



1 + (X(s, r 2 , n) + D(s, r 2 , r 1 ))" 1 X(s, r 2 , rf) 



[ e w — £ m) 



(95) 



and the total free energy of the tube is 

F C (R,5) =F£> W (R,5)-Fj°\5) 



fS ] w (s,R,S) = -k B T\og(det[B( S ,r 2 ,n)}) . 



(96) 



The required free energy Fc(R, S) is then given by Eqs. (73) and (74). Using Eqs. (52) 
and (46) it can be verified, as expected, that F^ W (R, 5 = 0) = 0. From Eq. (92) we 
find also that F^R, 8) =0 when 
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6. Evaluation of the Casimir energy 

In this section we evaluate the Casimir energy for the dielectric tube as a function of the 
inner radius R formed from a membrane of fixed thickness 5 and dielectric constant e^- 
The regions interior and exterior to the tube are water-filled with dielectric constant 
ew- The cross-section of the tube is shown in Figure 2 and the length of the tube is 
aligned along the z-axis. 

The result for the free energy of the tube, Fc(R, 5), is given in Eqs. (96) as a sum 
over the mode number n and an integral over the wave- vector, p, in the ^-direction. We 
evaluate the sum and integral numerically and present results in the next section for 
various values of 5. However, as is usual in many cases, the calculation is dominated by 
the Ultra- Violet (UV) properties of the integrand and an UV cut-off must be imposed 
to achieve a finite result. We examine the UV properties of the integral and calculate 
the leading divergent contributions analytically. These divergent contributions, which 
are regulated by the UV cut-off, agree with the prediction for them obtained from the 
full numerical calculation. We also verify that the R — > oo limit of F$ W (R,8) agrees 
with the result in the planar film case for a film of thickness 5. It is convenient to define 
the following constants which encode the dielectric properties of the system 

A = ■ , 7 = 1 = . (97) 

+ cm £m i — za 

After some algebra and use of the Wronskian identity Eq. (52) we find 

/ Inipr^Kniprf) I n (pri)K n (pr 2 ) \ 



(D(s,r 2 ,ri) + X(s,r 2 ,r 1 )y 



(98) 



y I n {pr 1 )K n (pr 2 ) I n (pr 2 )K n (pr 2 ) J 

Using the expression for the diagonal matrix X in Eq. (91) we find that B(s,r 2 ,ri), 
whose determinant is required for the evaluation of S^ ^ in Eq. (95), is 

/ 1 + 7^(^1)^(^1) -I n (pri)K' n (pr 2 ) \ 



B(s,r 2 ,rf 



(99) 



V I 'n{P r i)K n {w2) 1 - yI n (pr 2 )K' n (pr 2 ) J 

The important feature of B is that the on-diagonal elements are the separate 
contributions from the surfaces at r = r\ and r = r 2 and the off-diagonal terms are 
the contribution from the interaction between the surfaces. In particular, it will be 
shown later in this section that the off-diagonal elements fall off exponentially with the 
surface separation 5 like exp(— 2^J~(p 2 + m 2 )5). This fact has the consequences that the 
inter-surface interaction becomes negligible for large separations or for large wave- vector 
p. The corollary is that any Ultra- Violet divergences, which are due to large p behaviour 
of the integrand, arise solely from the separate surface contributions and that the inter- 
surface interaction gives a UV-finite contribution. We therefore explicitly separate the 
terms in F^] V (R,8) into these respective contributions. We have 

F ^ {R 1 S) = -g(A ri , A) + -g(Ar 2 , -A) + fc(n, r 2 , A, A) + m(A, A) , (100) 
Lk B T n r 2 
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where 

g(x,A) = f X dulog[l + A(I n (u)K n (u)y} , 

2tt - Jo 

h(r 1 ,r 2 ,A,A) = 

1 M dk ^ 4 A Vr i r 2 1^ (pn ) I n (pr x ) K' n (pr 2 ) K n (pr 2 ) 

2 Jo tv ^ ° g [ (l + Apr 1 (I n (pr 1 )K n (pr i y)((l-Apr 2 (In(pr2)K n (pr 2 ) / ) 

m(A,A) = -Jf-EMl-A 2 ). (101) 

2 Jo n V 

The contribution LksTg(Ar, A)/R is the free energy of an isolated cylinder of length 
L, radius i? and dielectric constant eu in a medium of dielectric constant 6w Thus 
the first two terms in Eq. (100) are the respective separate contributions of the inner 
and outer cylindrical regions that form the layer of thickness 5 = r 2 — r x ; the term 
LksTh(ri, r 2 , A, A) is the contribution from the interaction between the cylinders. As 
expected, the function g(x, A) diverges as x — > oo and so this term in the free energy 
must be regulated by taking a finite non-zero cut-off A = ir/a, where a is the UV 
cut-off length. Viewed as a Taylor expansion in A we find that the 0(A) term of g is 
independent of r and so in the free energy the contributions proportional to A cancel. 
This to be expected on physical grounds since by examining the limit of a diffuse system 
one can see that any term proportional to A must be a self energy term [25]. The term 
of order A 2 of g can be evaluated using Bessel function summation theorems [24] . This 
term is given by 

g 2 (x,A) = [ X duU 2 Y,iIn(u)Kn(u)f . (102) 

4"7T Jo „ 



We define R*(r 1 ,r 2 , <p) = \Jr\ + r\ — 2r\r 2 cos(0) and from [24] we have 

K Q {uR*{r x ,r 2 ,<f)) = Y, I n(ur l )K n {ur 2 )e m4> . (103) 



Then 



Thus 



K Q {2uR*{l,l,<t>)) =Y,WK n (u)e m4> =}► 

n 

2 sin(0/2) K' (2u sin(0/2)) = ]T [I n {u)K n {u)\ e m * . (104) 

n 

J2 [I n (u)K n (u)f = - / K\{2u sin(0/2)) sin 2 (0/2) , (105) 

n TX JO 

where the Bessel function identity K' Q {u) = —K 1 (u) has been used. By substitution into 
Eq. (102) and careful manipulation of the double integral we find that 

g{x,A) = --^-A 2 [61og(x) + 301og2 + 67-ll] + 0(A 4 ) + 0(l/x) . (106) 
zoo 

Similarly the finite contribution h(ri, r 2 , A, A) can be expanded and the 0(A 2 ) term 
evaluated. We have 

h 2 (n,r 2 ,A,A) = -A 2 [ A dpnr.Y P 2 In(pri)Upri)K(pr2)K n (pr 2 ) . (107) 

71 Jo „ 
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Using Eq. (103) we have 
Y^p 2 rMi)In{prf)K' n {pr 2 )K n {pr 2 ) 



n 



(108) 



Because the integral is convergent we may set A to oo. We find 
k 2 (r 1 ,r 2 ,A,A) 



= ^ir d uK l {uf r<v 

IX JO JO 




(109) 



with r\ — R — 8/2, r 2 = R + 5/2. After manipulation we find that 





(110) 



where y = 5/2R. This gives 



» a x a 2 3 a 2 r / 5 \ n , n 1 

h 2 (r 1 ,r 2 ,A,A) = — + log — +21og2-- 

v ' ' ' ; 25 2 64 R I B \2RJ B 2 



(111) 



6.1. R — > oo, 5 /ixeo? 

To calculate F C (R, S) in the grand canonical ensemble we must subtract from F^|y(i?, 5) 
the free energy F^, defined in Eq. (75), for a flat membrane of the same area and 
thickness. This has been calculated in previous work [4] but is it instructive to derive 
it directly from Eq. (100). In the limit R — > oo the arguments of all functions for p ^ 
become large and we find that the calculation is dominated by large n. The leading 
asymptotic results given in Eq. (A. 3) will be sufficient to compute F^ w (s,R,8) in the 
large R limit. 

From Eq. (A. 3) we have for large n that 



In the limit R — > oo it is better to define a new two-dimensional wave- vector k = (n/R,p) 
(i.e., ki = n/R, k 2 = p), since we then find that the R — > oo limit can be formulated in 
terms of functions with finite arguments. The measure is then 




1 



1 



(112) 




(113) 



Note that A = 2nLR is the area of the tube. 
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From Eq. (112) and using the definition of r)(z) in Eq. (A. 4) we obtain 



rj { ^(R + S/2))-rj{^(R 



8/2) 



The next correction is 0(5 3 /R 2 ) which is negligible in the R 

{i n {wi)K n {wi)y 

I'Swi)Uw2)K n ( Wl )K' n {w2) 



hSv'ih/h) = {kl + klf'H = k8 . (114) 

oo limit. Then we find 
0, (115) 
(116) 



(I n (pr 2 )K n (pr 2 ))' 

exp(—2k5) . 



4£> i r 1 r 2 

Using Eq. (115), we see immediately that the contributions from the individual surfaces 
vanish in this limit. The non-zero contributions then arise only from h(ri, r 2 , A, A) and 
trivially from m(A, A) in Eq. (101). On substitution in to Eq. (100) we find in the 
large R limit that 

d 2 k 

We consider this result in the limit 5 — > 0. The integral in Eq. (117) must be regulated 
with a UV cutoff k < A. The second logarithm in this equation can be expanded and 
the series in A 2 integrated term by term. If we assume that A5 3> then we find 



-\j^T2 [log(l"A 2 ) - log(l-A 2 exp(-2M 



(117) 



PF M = - - 



^ log (l - A 2 ) ^-L E 



,2m 



(118) 



25 2 — m- 

The result behaves like l/<5 2 but only as long as the assumption A5 3> holds since 
terms containing the factor exp(— A5) have been ignored. If all terms are kept then, of 
course, lim^ F<x>(8) = 0. 

On subtracting 2ti RLksTF^ from F^ w to obtain the grand free energy Fq we see 
that the first term in cancels the contribution from m(A, A) in Eq. (100) identically. 
We retain the second term in at 0(A 2 ) and it cancels a similar term in the evaluation 
of the integral for h 2 (ri, r 2 , A, A). This term is exhibited explicitly in Eq. (111). Putting 
our results together we find that k c defined in Eq. (4) is given by 



k c = 



A 2 
64 



3 log 




+ 61og2 + 3 7£ -4 



+ A 4 £?(A) , 



(119) 



where 7# is Euler's constant and the constant in the brackets is evaluated to be 
0.02954. . .. We note an important point which is that the log(i?) dependences from 
the functions g and h cancel exactly giving a leading order behaviours of F c ~ 1/R. 



7. Numerical Results 

In order to calculate the Casimir energy as a function of R and S we evaluate Fq(R, 5), 
defined in Eqs. (96) and (100), numerically. The free energy Fc{R,S) is normalised 
to zero for R = oo, and is defined in terms of the free energy F^ W (R, 5) which is 
normalised to be zero when the dielectric constant of the membrane, €m is set equal to 
that for water: €m = %■ Here, F^ W (R, 5) is given as a sum over n and integral over p 
of F$ w (s, R, 5) (s = (n,p)), itself defined in Eq. (96). 
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To carry out both the sum over n and the integral over p we use the VEGAS 
integration package [26] which is an efficient algorithm which uses importance sampling 
to do multidimensional integrals. Although we are dealing with a discrete sum over n 
it is easy to adapt the integrand so that it is a function of the continuous variable x 
through the relation 

n(x) = Rlnt(x) , x = x - 0.5(1 - sign(x)) . (120) 

Then n(x) takes integer values necessary for the summation whilst x is used as a 
continuous integration variable by VEGAS. Both the sum over n and the integral 
over p are done by efficient importance sampling techniques, and an accurate answer 
can be obtained. To impose the needed Ultra- Violet regulator or cut-off, we set 
k = (n(x)/R, p) and integrate over the region — n/a < k < —ir/a , k — |k|. 

The evaluation of the integrand poses some difficulties since, as we have seen in 
the previous section, the integrand is dominated by large values of k, and hence the 
arguments of the Bessel functions in the definition of F&(s, Eqs. (100) and 

(101), become very large indeed. In this case the function /„ (K n ) suffers from floating 
point overflow (underflow), which can be seen easily from the asymptotic forms given in 
Eq. (A. 3). However, in contrast the products over the Bessel functions which constitute 
each term in Eq. (101) do not suffer in this way. This also can be seen from Eq. (116) 
where the increasing and decreasing exponential behaviours of /„ and K n , respectively, 
compensate to give the behaviours exp(— 2k5). To construct a robust integrand we 
used routines for the full Bessel functions when k was sufficiently small and used the 
appropriate asymptotic form given in (A. 3) when either n or p, or both, became large. 
It was then possible to cancel the diverging and vanishing exponential factors against 
each other, so obtaining a well defined integrand computationally. 

We take ew/^M = 40 and evaluate Fc(R, 5) as a function of R in nanometers, and 
5 = 5(nm), 10 (nm) and for various values of the cutoff length a. Because there is no 
electrolyte the temperature dependence is purely in the factor of ksT multiplying our 
calculation. From Eq. (119) it is clear that, as is true in most applications, the Casimir 
energy is dominated by the UV cutoff behaviour and hence by the value chosen for a. It 
is not fully clear what the correct value for a should be since the microscopic properties of 
the membrane interface are not properly included in the analysis. Typically, we would 
expect a to be the scale of the inter-molecular spacing of the molecules forming the 
membrane or of water molecules. For this calculation a reasonable value is a ~ 0.5(nm). 
To test the validity of our UV analysis we first investigate how Fc(R, 6) behaves for 
very small a and we choose a = 0.05, 0.1, 0.2, 0.5(nm). The function B(A 2 ) receives 
contribution from the both the g and h terms in Eq. (27) with -B(O) 7^ 0. Note that 
there are no odd terms in A in the leading 1/R behaviour of Fc since to leading order 
one may set 5/R = (equivalently i?i = R 2 ) in the leading order behaviour of g and 
in the denominator of the second term in the logarithm of the integral defining h. This 
is a consistent parametrisation whilst 5 ^> a. The limit 5^0 must be taken carefully 
and when 5 < a the separation of F^ w in Eq. (100) into contributions from functions 
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78/82 


10 3 


-0.342 


-0.443 


0.123 


78/82 


10 2 


-0.244 


-0.346 


0.123 


0.6 


10 3 


-0.1361 


-0.1520 


0.123 


0.6 


10 2 


-0.0972 


-0.0162 


0.123 


0.2 


10 3 


-0.0151 


-0.0162 




0.6 


10 3 


-0.0038 


-0.0040 





Table 1. For various values of A and S/a we compare the prediction of Eq. (119) 
with the result of simulation and deduce a numerical value for B(A 2 ). Owing to 
small systematic errors in the numerical calculation of the Bessel functions there is a 
negligible discrepancy for very small A but £>(A 2 ) is seen to be a constant function 
from evaluations at larger A. We see that the result for Fc from Eq. (119) is in very 
good agreement with the full calculation. Various values of S and a were used but 
typically (5=1 — 10 (nm) 

g and h is not useful since h develops the compensating UV divergence to that in g and 
we find lim^o F c = 0, as expected; in essence, the larger of (S : a) acts as the UV cut-off 
on the integral defining h. 

In Table 1, for various values of A and S/a, we compare the prediction of Eq. (119) 
with the result of numerical evaluation and deduce a numerical value for B(A 2 ). Owing 
to small systematic errors in the numerical calculation of the Bessel functions there is 
a tiny discrepancy for very small A but B(A 2 ) is seen to be a constant function from 
evaluations at larger A and we see that B(0) is plausibly 1/8. 

The physical value of the UV cut-off length can only be determined 
phenomenologically in this model. This is because the model is an effective field theory 
in which the dynamics of the molecular electric dipoles is described by the dielectric 
constant which is a static long-range parameter. The field modes with large-/c and n 
probe the static short distance properties of the model and so a more refined field theory 
is needed for these scales. It is unclear whether the molecular nature of the lipid has an 
effect on the UV cutoff but it would seem most likely that the effective value of ew at 
short scales are dominant in this calculation. The effect is encoded in the value of a. 

8. Conclusion 

In this paper we have developed the theory for a general approach to the calculation 
of the electrostatic free energy for a system of symmetrically layered electrolytic 
membranes. The definition of symmetrical is that the Laplacian is separable and that 
the coordinate direction normal to the layers can be interpreted in terms of a Euclidean 
time variable, t. The x partition function for a layer bounded by surfaces defined by 
t — ti and t = tj, ti < tj can then be written in terms of the Feynman (Euclidean) 
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time evolution kernel from ti to tj by invoking the well-known connection between this 
formalism and statistical mechanics. The method is a general extension of our work 
concerning flat membranes [4], [27], [22], [8] and allows the full interacting Sine-Gordon 
field theory in Eq. (7) for electrolytic layers to be studied systematically, including the 
perturbation theory in the coupling constant g defined by g = Ib/Id, where Is and Id 
are the Bjerrum and Debye lengths, respectively. Geometries of interest to which this 
approach applies include cylindrical and spherical ones. In this paper the general theory 
developed in section 3 is applied to the system of cylindrical layers each filled with a pure 
dielectric medium where the dielectric constant differs between layers. The analysis of 
this system is based on the free harmonic field theory which is exactly soluble and which 
is the theory about which a perturbative expansion for the effect of fluctuations takes 
place. In the succeeding sections the particular problem of the free energy of a tube 
of dielectric material submersed in water is studied and the Casimir force calculated. 
The tube is a simple model for the t-tubule formed from lipid membrane in particular 
kinds of muscle cell and the object of the calculation is to investigate the size and 
form of the Casimir force and whether it can act to stabilise the tube against bending 
stress that acts to decrease the tube curvature and hence increase its radius. In order 
to calculate the relevant free energy we assumed that the tube was connected to an 
infinite reservoir of flat membrane and so worked with the grand ensemble describing 
this system. The assumption about which ensemble is the relevant one is crucial to 
answering questions concerning tubule stability since in the grand ensemble the nature 
of the lipid or membrane reservoir determines its bulk free energy and hence affects the 
energy of conformation of the tube; a flat membrane reservoir will differ for a reservoir 
which stores spare lipid at the end of the tubule essentially as crumpled membrane. 
Clearly also, if there is no reservoir at all, so that the ensemble is canonical, then any 
increase in the surface area of the tube as it expands will lead to a lower density of lipids 
in the surface and a concomitant increase in the surface energy over and above that due 
to the bending stress; this is a fundamentally different situation to the one we assume. 

In the case studied here we see from Table 1 that for a lipid bilayer tube in water 
with 5 = 10(nm) and a = O.l(nm) we find kq = 0.346. It is possible to include 
the contributions from the modes with non-zero Matsubara frequencies, a calculation 
in progress, but we can expected at most a factor of two or so enhancement on past 
experience of similar calculations ([2]), and so k c ~ 1 is a likely largest value . These are 
at the lowest end of values for Kb for known lipid bilayers in water ([10, 11]). However, 
Wiirger ([11]) calculates k b for surfactant films, analysing the role of hydrophobic tails, 
as a function of the tail length and the area per molecule, and finds a wide range 
of values for k b (k_b = 7T«, with k from ref. ([11])) including values small enough, 
corresponding to soft interfaces, to accommodate our result. Thus it is conceivable that 
there can be small tubes formed from soft membranes in water for which the bending 
forces tending to expand the radius are compensated by the Casimir attraction, and the 
tube is stabilised by sub-leading 0(1/R 2 ) forces. 
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Appendix A. Asymptotic behaviour of modified Bessel functions 

For completeness we include the asymptotic behaviours of the modified Bessel functions 
I n and K n as given in ([28]). For z — > 0: 

T (A {z/2r 

K n (z) ~ \n\{z/2)- n (n>0) 

K (z) ~ -{log(z/2)+ 7 }J (*). (A.l) 
For z — > oo: 



V ; V27TZ \ 8z 2!(8z) 2 J 

( (//-!)(// -9) 1 



where /x = 4n 2 . 
For n — > oo: 



In{nz) ^7S(TT^l 1 + £^) 

V2^ z 1 fc tt n k J 



where 



f = 7rr?' " = ^ +l0g TT7rr7' (A ' 4) 

and Mfc(t), Vk(t) are polynomials of order 3k and are even (odd) if k is even (odd) [28]. 
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